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Abstract 

We study the effect of random dilution of attraction centers on electron pairing 
in a generalization of the negative U Hubbard model where the interaction C/j, 
defined on the two dimensional lattice, can be — 1?7| and with probability c and 
1 — c depending on the lattice site i. Using the determinant quantum Monte Carlo 
approach, we find the critical concentration cq of negative centers which leads to 
superconductivity, and show how the evolution of the local density on the two types 
of sites lends insight into the formation of global phase coherence. 
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1 Introduction 



The behavior of High Temperature Superconductors (HTS) is characterized by 
a strong dependence on doping. For example, changing the Sr concentration in 
La2_xSr2,Cu04, or the O concentration in YBa2Cu306+5 tunes these materials 
from antiferromagnetic insulators, to spin glasses, and finally to superconduc- 
tors. The short coherence length and relatively low superfiuid density in the 
underdoped region emphasize the existence of phase fluctuations and enhance 
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the competition between superconductivity and other long range order phases 
[1,2]. In such systems, where different ground states emerge as the particle den- 
sity is tuned, it is often important to explore the possibility inhomogeneous 
mixtures of phases; for example, 'antiferromagnetic' and 'superconducting', 
at commensurate and incommensurate filling. Indeed, there has already been 
considerable study of whether the doped holes are uniformly distributed or 
whether they accumulate preferentially in two-dimensional phase separated 
regions[3] or in one- dimensional stripes[4]. 

There have been many direct studies of doping dependent transitions in corre- 
lated electron Hamiltonians like the repulsive Hubbard model, which provide 
microscopic models of IITS[5,4]. These are numerically challenging since they 
attempt to understand the origin of HTS at the same time as addressing ad- 
ditional subtle issues like inhomogeneous phases. A more simple, if somewhat 
less fundamental, way in which related questions can be addressed is to con- 
sider models, like the —U Hubbard Hamiltonian, which have an attractive 
interaction 'built in'. This separates out the problem of how superconducting 
pairs form in disordered, short coherence length systems in which the particle 
density is fixed but the interactions are diluted, from the complexities of the 
origins of the attractive interaction. The key underlying question is how super- 
conductivity spreads through a system as the local attraction between charge 
carriers is strengthened. Indeed, such centers of attraction played important 
role in the several models of HTS [6-14]. Another possible experimental re- 
alization is the percolative, granular-like transition to superconductivity via 
hydrogen doping which has been reported in Eui.5Ceo.5NbSr2Cu20io[15], a 
ceramic with the superconducting critical temperature Tc— 29 K. 

In this paper, wc shall adopt this approach and study a simple model which 
allows us to address the issue of superconducting transitions tuned by changing 
the concentration of attraction centers. Specifically, we study a two-dimensional 
negative U Hubbard model which has an on-site attraction on some fraction 
of its lattice sites and ask what fraction of sites is needed to be attractive for 
the system to supcrconduct. Recently, this issue has been addressed by Litak 
and Gyorffy [16] who examined the problem of superconducting percolation 
by using Hartree-Fock-Gorkov decoupling as an approximation for interacting 
electrons and Coherent Potential Approximation (CPA) to treat randomly 
distributed centers. They have found the critical concentration of negative 
centers cq = n/2, where n denotes a band filling in the limit of large \U\ inter- 
action. For a small concentration of centers c < Cq the system becomes normal 
because every site, where —\U\ interaction is present, is doubly occupied. 

The effect of dilution of the interactions in the attractive U Hubbard model at 
half-filling, n = 1 can also be inferred from the Quantum Monte Carlo work 
of Ulmke et a/. [17], who investigated the suppression of antiferromagnetism 
hy U — impurities in the repulsive Hubbard model ( C/ > 0) at n = 1. The 
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critical concentration fc of impurities which destroy the antiferromagnetism 
was estimated as ^ 0.45. Since a particle-hole transformation maps the 
attractive and repulsive U Hubbard models onto one another, at half-filling, 
this work immediately identifies the critical concentration cq for suppression 
of superconductivity in the attractive Hubbard model, cq — 1 — fc — 0.55. 

Further analytic work on the disordered attractive Hubbard model, for ex- 
ample based on numerical solution of the Bogliubov-de Gennes equations, is 
contained in [18-20]. However, there have to date been no Quantum Monte 
Carlo (QMC) studies for the doped case where the density is not precisely 
half-filled and particle-hole symmetry is broken. In this paper we examine 
the problem of percolating superconductivity using QMC in the negative U 
Hamiltonian changing the concentration of negative centers c as well as band 
filling n[16]. The use of QMC enables us to go beyond the hmitations of mean 
field approximations and investigate the effect of charge and pairing potential 
fiuctuations on an equal footing. 

The paper is organized as follows. After an introduction in Sec. 1 we present 
the microscopic model and formulate the problem (Sec. 2). Section 3 defines 
the QMC Monte Carlo method used in siniTiIations and shows the key numer- 
ical results: the prediction of superconductivity percolation and the critical 
concentration cq. The paper ends with conclusions. 



2 The microscopic model 



We consider the random negative U Hubbard Model [16] defined by the fol- 
lowing Hamiltonian: 

H — — tijC^Cj^ -\- — Ui C^^Cif^cl_^Ci^fj — /X C^a^io- ■ (l) 
ija ia ia 



Here i and j label sites on a two-dimensional square lattice, tij is the hop- 
ping integral connecting only nearest neighbour i and j's, fi is the electronic 
chemical potential, Cj^, Ci^ create and annihilate, respectively, electrons at the 
single site i with spin a, and the coupling constant. 



— \U\ with probality c 

with probability 1 — c 



(2) 



We are interested in determining, away from particle-hole symmetry, the criti- 
cal concentration cq for the U centers such that for c < cq the configurationally 
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Fig. 1. The two dimensional lattice with random electron interaction U = or — |C/|. 



averaged, superconducting long range order parameter vanishes. At the same 
time, we would like to obtain some insight into the mechanism of the formation 
of long range order as one moves to c > Cq. 

Such questions can be addressed by computing expectation values of operators 
like the Hamiltonian, the kinetic energy and double occupation, and pair-pair 
correlation functions. Such expectation values take the general form. 



and can be evaluated by a standard procedure [21-24] which involves writing 
a path integral expression for the imaginary time evolution operator e~^^ and 
performing a stochastic integration over a set of auxiliary fields introduced to 
treat the electron-electron interaction exactly. This will be described in greater 
detail in the following section. Typically, we can study system of a few hundred 
electrons, far larger than the competing exact diagonalization approach, but 
still small enough to require careful treatment of finite size effects, as we shall 
discuss below. 

We will focus especially on the equal time pair correlation function, which is 
defined by 

Ps{j-l)^<^At >: where A+ = c+ c+ (4) 
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Here creates a pair of electrons at lattice site i. The long range order in 
Ps{l) can be probed by calculating its structure factor, 



A finite positive value Sp > 0, in the limit of large lattice size L. signals 
superconductivity, otherwise the system is normal. For our disordered system 
(Eq. 2, Fig. 1) we have to do a disorder average over many realizations of 
random configurations of the attractive U centers. 



3 Details of Quantum Monte Carlo Approach 



In order to describe the quantum Monte Carlo method [21] in more detail, let 
us first rewrite the Hubbard Hamiltonian in a way which puts the interaction 
term in a particle-hole symmetric form. 

ija ia i 



The trace over the fermion degrees of freedom cannot be performed analyt- 
ically due to the quadratic interaction of term U. To reduce the problem to 
a quadratic Hamiltonian we discretize the imaginary time (3 = LAt (Eq. 3) 
and employ the Trotter approximation [25,26] to decompose the full imaginary 
time evolution. The partition function can be expressed as. 





L 











(7) 



where K includes the quadratic part of H (Eq. 6) 



ija ia 



(8) 



while P the on-site interaction part 

P = Ef^.(^.T-^)K-^) ■ (9) 
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Then wc apply the discrete Hubbard-Stratanovich transformation to decouple 
the attractive interaction Ui < [22-24] 



2 Si=±l 



where 

cosh (ArA) = e^'l^l/^ (11) 



The replacement of the quartic interaction term between up and down spin 
electron densities by an expression in which the densities couple to a Hubbard- 
Stratonovich field makes all terms in the exponential of Eq. 7 quadratic in the 
fermion operators. The trace can then be performed analytically, leaving a 
product of determinants which depend on the Hubbard- Stratonovich field. 
Because the up and down spin species couple to the field with the same sign, 
these two determinants are identical, and hence their product is always posi- 
tive. So, unlike the case of the repulsive model where — riii couples to the 
field, there is no sign problem in the attractive model, even in the presence of 
disorder in the interaction. 
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Fig. 2. Finite lattice scaling of the structure factor S^. Each point represents the 
average for 20 different disordered lattice configurations (for 4 x 4, 6 x 6 and 8x8) 
with relatively small error. 



4 Percolation of Superconductivity 



Using the above procedure we have simulated the Hamiltonian (Eq. 6) and 
calculated the pair structure factor Sp (Eq. 5) as a function of lattice size L. 
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Spin-wave theory suggests that in the ordered phase the finite size correction 
for Sp should be proportional to the inverse linear system size 1/L[27]. Figure 
2(a) presents the results of finite size scaling of Sp for various concentrations 
of U sites c. The calculations were done for U — —At and P — 1/T — 8/t, 
and the average number of electrons per site was chosen to be n = 1. One can 






Fig. 3. Charge on U site nu and the structure factor of pair correlation function Sj 
versus system filling n for a 8 x 8 lattice. 
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easily see that there exists a critical concentration cq ~ 0.75. 



To see how Cq depends on n in Fig. 2(b) the concentration of U centers is 
kept fixed at c = 0.750, and we have plotted Sp for various system fillings 
(n = 0.75, 1.00 and 1.25). Here n = 1 corresponds to the case at the percola- 
tion threshold (c = Cq). It is clear that for a larger average charge n = 1.25 
the system is not superconducting. A different situation emerges for n = 0.75 
where Sp{L — > oo) appears to be slightly above 0, that is, the system is super- 
conducting. Clearly, disorder in the interaction U has broken the particle-hole 
symmetry. In fact, this is connected with the modification of the bare chemical 
potential at U sites, in the quadratic part of Hamiltonian K (Eq. 6,8) 



/.-t//2 for a. = -\U\ ^^^^ 
II for Ui — 



In some previous work on the repulsive model[17, 22-24], disorder was inserted 
in a way which preserved particle-hole symmetry by using a separate shift in 
chemical potential depending on the interaction strength on each site. This 
enabled simulations at lower temperature since it eliminated the sign problem, 
but also was somewhat less physically motivated. 

Interestingly, the structure factor associated with the pairing correlation func- 
tion Sp{L — >• oo) for c = 1 is smaller than the one for c = 0.875. This is due 
to the existence of a competing charge density wave (CDW) phase at n = 1. 
A similar effect was found for the system with disordered site energies and 
the same attraction U on every lattice site [22,28]. Thus, in addition to phase 
fluctuation of the pairing potential[18,19], it is clear that the distribution of 
charges on the lattice (charge fluctuations) is also important for the destruc- 
tion of superconductivity. Indeed, previous work suggests that percolation of 
the occupied sites is a crucial phenomenon in this process [16]. 

Examining this effect further, one can look for the change of the pairing factor 
correlation function Sp and simultaneously the average occupation of U sites 
(nu) (where U — —\U\). Figures 3a-c show these quantities as a function of 
average site occupation n for a flnite system 8 x 8 for c = 0.875 (a), c = 0.625 
(b) and c=0.375 (c). respectively. Note that Sp docs not display particle- 
hole symmetry. The dotted lines in Fig. 3(a-c) represent the occupation of 
attractive sites for a simplified Hamiltonian, where K reads 

K ^-tY, ctcja - Yl /^^i^ ■ (13) 

ija ia 



It is easy to see that in case of the full Hamiltonian (Eqs. 2,6) the system 
displays charge fiuctuation {nu 7^ n) and comparing nu curves to the dotted 
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lines one can see that charge fluctuations are stronger for smaller c. In our 
case, for large enough interaction U, the small concentration of U centers 
leads to double occupation in almost all attractive sites. This makes it harder 
for electrons to hop to another [/ 7^ site. Nevertheless the system can be 
metaUic because the other electrons which occupy the U — site can easily 
hop to the other U = sites. This result is consistent with the results obtained 
by Litak and Gyorffy [16]. The small dip visible in the pairing structure factor 
Sp around n = 1 (in Fig. 3a) can be associated with the creation of a CDW 
order. This effect is not present in Figs. 3b-c as CDW expected to disappear 
faster than superconductivity [22,28]. 



5 Conclusions 

We have examined the problem of percolating superconductivity in the con- 
text of a random U Hubbard model. We have studied the case where Ui is 
— \U\ and with probability c and 1 — c respectively by the quantum Monte 
Carlo method. For the half-filled system n = 1 we obtain the critical concen- 
tration Co = 0.75. As we expected the transition is due to phase and charge 
fluctuations and depends strongly on the density n. In that sense our results 
agree qualitatively with those obtained by Litak and Gyorffy [16]. Moreover, 
the present approach is exact as opposed to earlier work [16] employing the 
mean fleld approximation. We flnd that the onset of superconductivity in this 
diluted, attractive Hubbard model is linked to percolation of the order param- 
eter and charge fluctuations. 
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